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Abstract 

We compute fragmentation corrections to hadroproduction of the quarkonium states J/i), Xcj, 
and V’(2<S') at leading power in where me is the charm-quark mass and px is the quarkonium 

transverse momentum. The computation is carried out in the framework of nonrelativistic QCD. 
We include corrections to the parton-production cross sections through next-to-leading order in the 
strong coupling ag and corrections to the fragmentation functions through second order in ag- We 
also sum leading logarithms of p‘^/m‘1 to all orders in perturbation theory. We find that, when we 
combine these leading-power fragmentation corrections with fixed-order calculations through next- 
to-leading order in ag, we are able to obtain good fits for px > 10 GeV to hadroproduction cross 
sections that were measured at the Tevatron and the LHC. Using values for the nonperturbative 
long-distance matrix elements that we extract from the cross-section fits, we make predictions for 
the polarizations of the quarkonium states. We obtain good agreement with measurements of the 
polarizations, with the exception of the CDF Run II measurement of the prompt polarization, 
for which the agreement is only fair. In the predictions for the prompt-J/V’ cross sections and 
polarizations, we take into account feeddown from the Xcj and '4>{2S) states. 

PACS numbers: 14.40.Pq,13.88.-Fe,13.87.Fh,12.38.Bx 
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I. INTRODUCTION 


In recent years, corrections to inclusive quarkonium production cross sections and polar¬ 
izations through next-to-leading order (NLO) in the strong coupling have been computed 
for both hadroproduction [1-6] and photoproduction [7-9]. These computations have been 
carried out in the context of the nonrelativistic QCD (NRQCD) factorization conjecture [10], 
which states that the inclusive production cross section to produce a quarkonium 77 in a 
collision of particles A and B can be written as 

dcTA+B^B+X = E d0'A+B^QQ(n)+x(^^('n)}- ( 1 ) 

n 

Here, the da A+B^QQ{n)+x fhe short-distance coefficients (SDCs), which can be computed 

in perturbation theory and which correspond to the production of a heavy quark-antiquark 
pair QQ{n) in a specihc color and angular-momentum state n. The {0^{n)) are NRQCD 
long-distance matrix elements (LDMEs), which parametrize the nonperturbative part of the 
production process. 

Because the LDMEs have a known scaling with v, the heavy-quark velocity in the quarko¬ 
nium rest frame [10], the sum in Eq. (1) can be regarded as an expansion in the small 
parameter v. ~ 0.3 for the J/ip.) In present-day phenomenology, the sum in Eq. (1) 
is truncated at relative order v^. For H = J/ip or H = ■^(25'), the truncated sum involves 
four LDMEs: and where the expressions 

in parentheses give the color state of the QQ pair (singlet or octet) and spin and orbital an¬ 
gular momentum in spectroscopic notation. Here, ip stands for J/ip oi ip{2S). For H = Xcj, 
the truncated sum involves two LDMEs: and where the LDMEs 

for the Xci and Xc 2 states can be related to the LDMEs for the Xco state by making use of 
the heavy-quark spin symmetry [10], which is valid up to corrections of relative order 

Since the color-singlet LDME for quarkonium production is related to the 

color-singlet LDME for quarkonium decay, it can be determined in lattice QCD, from po¬ 
tential models, or from the ip decay rates into lepton pairs. On the other hand, it is not 
known how to compute the color-octet production LDMEs from first principles, and they 
are usually hxed by comparisons of NRQCD factorization predictions with measured cross 
sections. 

Even at the level of NLO accuracy in the theoretical predictions, it is not possible to 
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achieve a fully consistent description of the existing J/'ip production data within the NRQCD 
framework. For example, one can £t the hadropro duct ion cross-section data [11, 12] and 
polarization data [13-15] simultaneously [4], but the LDMEs that are obtained yield a pre¬ 
diction for the photoproduction cross section that is larger than the HERA data from the HI 
Collaboration [16, 17] by factors of 4-8 at the highest value of pt at which the cross section 
has been measured [18]. On the other hand, one can fit the predictions for the hadropro- 
duction and photoproduction cross sections to the experimental data [5], but the LDMEs 
that are obtained lead to predictions of large transverse polarization in hadroproduction 
at large pt, in disagreement with the experimental data [5]. In addition, it was fonnd in 
Ref. [19] that the pc production data that were measnred by the LHCb Collaboration [20] 
are incompatible with the LDMEs that were extracted in Ref. [5] from hadroproduction and 
photoproduction cross-section data. Although one can describe the pc production data by 
using the LDMEs that were extracted in Ref. [3], there is a very large cancellation between 
the contributions from the and channels [21, 22], and, hence, the remainder may 
be strongly dependent on uncertainties from nncalculated higher-order contribntions. 

These difficulties provide motivation for calculations of quarkonium production cross 
sections beyond NLO accnracy in An approach that simplifies compntations beyond 
NLO in as is to compute rates at leading power (LP) or next-to-leading power (NLP) in 
where rric is the charm-quark mass and Pt is the qnarkoninm transverse momentum. 
LP contributions can be factorized into semi-inclusive partonic cross sections to prodnce a 
specific single parton convolved with one-parton fragmentation functions (EEs) [23]. NLP 
contributions can be factorized into semi-inclusive partonic cross sections to prodnce two 
specihc partons convolved with two-parton FFs [24]. Calculations of these fragmentation 
contributions, at any given order in ag, are mnch simpler than a full fixed-order calculation. 
Fnrthermore, the LP- and NLP-factorization frameworks are natural ones within which 
to resum large logarithms of Of course, because the LP and NLP contributions 

represent the leading and first snbleading terms in an expansion in powers of ml/p'^, one 
would not expect them to be valid unless pr is significantly greater than uic- 

In Ref. [25] it was fonnd that LP contribntions beyond NLO in a* are important in J/ip 
hadroproduction. With the inclusion of these contributions, the LDMEs that are extracted 
from the prompt hadroproduction cross sections alone yield predictions for the J/ip polar¬ 
ization at large pt that are near zero and are in agreement with the experimental data [25]. 
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One deficiency in the analysis of Ref. [25] is that it does not take into account the effects of 
feeddown from the Xcj and states to the J/t/’- 

In this paper, we remedy that dehciency and extend the application of the LP- 
factorization approach by computing LP-fragmentation contributions to direct J/t/’, Xcj, 
and production. We extract LDMEs by htting to the Tevatron and LHC production 

cross sections, and we use those LDMEs to predict the J/t/’, Xcj, and polarizations. 

Our predictions for the prompt J/ijj and polarizations agree well with the existing 

high-pr LHC data, but the prompt J/ijj polarization is in only fair agreement with the high- 
Pt Tevatron Run H data. Our predictions for the Xcj polarizations will be tested soon at the 
LHC. While the results in this paper do not resolve the discrepancies between the NRQCD 
predictions and the J/ip photoproduction and pc hadropro duct ion data, they do provide a 
consistent description of the existing spin-triplet charmonium hadropro duct ion data at high 
Pt- 

The remainder of this paper is organized as follows. In Sec. H, we discuss the form of the 
LP corrections that we compute. Section HI contains the details of the calculation of the 
LP SDCs. We combine the LP and NLO results for the SDCs in Sec. IV. In Sec. V, we £t 
our predictions for the hadropro duct ion cross sections to the data, obtaining values for the 
LDMEs. We use these values for the LDMEs to make predictions for cross-section ratios 
and polarizations in Sec. VI. Finally, in Sec. VH, we summarize and discuss our results. 


II. CORRECTIONS TO QUARKONIUM PRODUCTION AT LEADING POWER 
IN pj ^ 


The contribution of leading power in pt to a quarkonium production cross section is given 
by the LP-factorization formula [23] 


A+B^QQirp+X^P) ~ 


j dz da 


A+B^i+x{Pi — P/z, Pf). 


( 2 ) 


Here, d&A+B^i+x is the semi-inclusive parton-production cross section (PPCS) for hadrons 
A and B to produce parton i, and D^^QQi^n) is the FF for parton i to fragment into the QQ 
pair with quantum numbers n. p is the momentum of the QQ pair, which is taken to be 
lightlike by neglecting the heavy-quark mass, and p* is the momentum of parton i, which is 
taken to be lightlike by neglecting the parton mass, pf is the factorization scale. 
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As we will describe in more detail in Sec. Ill, the PPCSs and the FFs have been calculated 
to order and respectively. Hence, we write them as 


do'A+B^i+x — + 0{ag), 

Di^QQin) — + 0{ag). 


(3a) 

(3b) 


As we have already mentioned, the SDCs for both unpolarized and polarized quarkonium 
production have been computed through NLO in Ug, which is order a^. In this paper, we 
extend these order-af calculations by combining existing calculations of the PPCSs through 
order and existing calculations of the FFs through order to obtain a partial calcu¬ 
lation of the order-a^ (NNLO) contributions to the LP SDCs. Furthermore, we calculate 
corrections to the LP SDCs involving leading logarithms of p'^/ml to all orders in by 
solving the Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) evolution equation [26-29]. 
Because this calculation of the LP SDCs accounts only partially for corrections of order af, 
we expect uncertainties from uncalculated corrections to be of order a^. However, these un¬ 
calculated corrections will not contain any enhancements from leading logarithms of Pr/ml. 

Part of the LP-fragmentation contribution through order is already included in the 
NLO SDCs, namely. 


c^^NLo(p) = ^ = viz, fif) 


+ I dzJ 2 »t {pi = pjz, {z, pf) 


( 2 ) 


= P/z, Pf)Dlll ^ {z, pf) 


( 1 ) 


i^QQ{n )' 


( 4 ) 


Hence, when we combine the SDCs through NLO in ag and the LP-fragmentation contri¬ 
butions, we must subtract the contributions in Eq. (4) in order to avoid double counting. 
Following Ref. [25], we compute 

dgLP+NLO ^ da^i^p 

dpT dpT dpT dpT ’ 

where da^i^o/dpx is the SDC through NLO in ag. The expression (5) takes into account, 
without double counting, the complete calculations through NLO in ag and also the addi¬ 
tional LP corrections beyond NLO that we have mentioned. 
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III. COMPUTATION OF THE LP SHORT-DISTANCE COEFFICIENTS 

In this section we describe the details of the computation of the PPCSs and FFs that 
enter into the LP short-distance coefficients in the LP factorization formula (2). 

We take nic = 1.5 GeV. We use the CTEQ6M parton distribution functions and the 
two-loop expression for with n/ = 5 quark flavors and = 226 MeV. We set the 

renormalization scale and the factorization scale fif for the both parton distribution 
functions and the FFs to be ttit = \/Pt 4m^. In order to resum leading logarithms of 
Pt/^1, we evolve the FFs from the scale po = 2mc to the scale pf = mr ~ Pt- We take 
the NRQCD factorization scale to be p\ = rric- In the calculation of the PPCSs and the 
evolution of the FFs, we take n/ = 3 active-quark flavors. That is, we ignore contributions 
from virtual or initial heavy quarks. 

A. Parton production cross sections 

The PPCSs through order were computed in the modihed minimal-subtraction (MS) 
scheme in Refs. [30, 31]. We carry out numerical computations of the PPCSs through order 

by making use of the computer code that was written by the authors of Ref. [30]. 

The PPCSs are computed as a function of pr, y, and 2 ; = p'^/pt = Pr/PiT, where pt is the 
transverse momentum of the QQ pair, y is the rapidity of the QQ pair in the hadron center- 
of-momentum frame, and piT is the transverse momentum of the specihc parton that is 
produced in the semi-inclusive partonic scattering process. Here, we have written in terms 
of the transverse momenta by using the fact that, in the LP approximation, one can ignore 
the invariant mass of the QQ pair. The maximum value of PiT is kinematically constrained, 
and, so, the PPCSs vanish for z < zo = + e“^), where i/i is the center-of-mass 

energy. 

B. Fragmentation functions 

In this paper we take into account FFs through order a^, which are available for frag¬ 
mentation of both gluons and quarks into polarized and unpolarized QQ pairs. A summary 
of FFs that we use in our calculation can be found in Ref. [32] and Ref. [33] for unpolarized 


7 



and polarized QQ pairs, respectively. We give a detailed description below of the sources of 
these FFs. 

The gluon FF Dg^qq^^n) for n = was calculated for both unpolarized and polarized 
hnal states at order (LO) in Ref. [34] and at order (NLO) in Refs. [32, 35]. The 
gluon FF for n = was calculated at order (LO) in Refs. [36, 37]. The gluon FFs 
for n = were calculated at order (LO) in Refs. [34, 37] for unpolarized hnal states 
and in Ref. [33] for polarized hnal states. The gluon FFs for n = ^Pj^^ were calculated at 
order (LO) in Ref. [34] for unpolarized hnal states and in Refs. [33, 38] for polarized hnal 
states. 

The situation for quark FFs Dg^qq^^^ with n an S'-wave state is rather complicated, as 
there are several independent calculations, some of which do not agree. Let us distinguish 
three cases: (i) g 7 ^ Q, in which case, n = (ii) q = Q and n = (hi) q = Q and 
n = The quark FF for case (i) for an unpolarized hnal state was calculated at order 
(LO) in Refs. [32, 39, 40], whose results all agree. The quark FF for case (i) for a polarized 
hnal state was calculated at order (LO) in Refs. [33, 39, 41]. The results in Refs. [33, 39] 
agree with each other, but disagree with the result in Ref. [41]. The results in Refs. [33, 39] 
have since been conhrmed by the author of Ref. [41]. The quark FF for case (ii) for an 
unpolarized hnal state was calculated at order (LO) in Refs. [32, 39, 40, 42]. The results 
in Refs. [32, 39] agree with each other and disagree with the results in Refs. [40, 42]. We 
use the results in Refs. [32, 39] in this paper. The quark FF for case (ii) for a polarized 
hnal state was calculated at order (LO) in Refs. [33, 39, 41], whose results agree. The 
quark FF for case (hi) for an unpolarized hnal state was calculated at order (LO) in 
Refs. [32, 39, 43], whose results agree. The quark FF for case (hi) for a polarized hnal state 
was calculated at order (LO) in Refs. [33, 39, 41], whose results agree. 

The quark FFs for n = ^Pj^' and n = were calculated for the unpolarized 

and polarized cases at order oi\ (LO) in Ref. [40]. 

The gluon FF Pg^qqqg\A^^ was calculated at order (LO) in Refs. [44, 45]. Because the 
contributions to the FF in the channel begin at order a^, we do not include them in 
our LP-fragmentation calculations. However, we do use the LO FF for the channel to 
estimate the size of the uncalculated LP-fragmentation contributions for that channel. 



C. DGLAP equation 


At leading order in as, the DGLAP equation is given by [26-29] 


d (Ds{fif)\ _ asii^f) ( Pqq 


gq 


dlogfif 27 r 


Pn. 


Pr, 


Pgil^f) 


( 6 ) 


99 99 

where Dg = Ds = E/ + L'g^^QQ(n)], / is the light-quark or light- 

antiquark flavor, the Pij are the splitting functions for the FFs, and nj is the number of 
active light-quark flavors. The symbol 0 represents the convolution 


{f<^ 9 ){z)= dx dyf{x)g{y)5{xy-z) = 


dx 


X 


f{z/x)g{x) = 


dx 


X 


The splitting functions are given by 
Pggiz) = 2 Ca 


+ + 4^ -z) + ^5(1 - z) 


f(x)g{z/x). (7) 


(8a) 


Pgq{.z) — Cp 


i + {i-zy 


Pggiz) = TP[Z^ + {I - Zf], 


P„(z) - Cr 


1 + ^" ^ 


(8b) 

(8c) 

(8d) 


where 


Ga 

Cf 

Tf 

bo 


= N,, 


E -1 

2N, ’ 

1 

2 ’ 

11 2 


= —Nr, 


;nf, 


(9a) 

(9b) 

(9c) 

(9d) 


and Ac = 3 is the number of colors. 

As is well known, an analytic solution to Eq. (6) can be obtained in Mellin space. The 
Mellin transform of a function / is defined by 

/(A) = iMf)iN) = [\zz^-^fiz), (10) 

Jo 

where we use a tilde ( ) to denote objects in Mellin space. The Mellin transform of the 
convolution in Eq. (6) is an ordinary product: 


lM{f ® g)](]V) = {Mf)(N) X {Mg)(N). (11) 
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Hence, Eq. (6) can be diagonalized by taking the Mellin transform. Using the one-loop 
evolution of as 


d 


_ _ 2 / \ 
dlog/ij dTT ^ dasifif)^ 

one obtains the following solution of the DGLAP equation: 

2P+/bo 


d 


Ds{N, fXf) 
Dg{N, /i/) 


M: 


(ho) 

as{i^f) 


+ M_ 


/ \ \ 2P /bo 

Q^s(/Xo) ^ 

as{fif) 


where 


M± = ±- 


Pgg-P^ 2nfP 


p+ _ p- 


P 


<13 


f^gq 

P — P^ 

^99 ^ 


( 12 ) 


Ds{N, ho) 

Dg{N , Po) 


(13) 


(14) 


and 


P^ = - 


^99 + Pqq V (^99 ^ 99 ^ + ^^fPqgPgq 


(15) 


The evolved FFs in ; 2 -space can be obtained from Eq. (13) by applying the inverse Mellin 
transform 


^C +200 

Di^QQin){z,^if) = — j dN z-^Di^QQ^^^{N,fXf), (16) 

ZTTt J c—ioo 

where the real number c is chosen so that the integral over N follows a contour that lies to 
the right of all of the poles of jj,f). 

We resum the leading logarithms of by choosing the evolution scales po = 2mc and 

f9f = rriT ~ Pt- In this paper, we compute the integral over N numerically, using analytic 
expressions for the Mellin transforms of the FFs at the scale po = 2mc- 

There is a difficulty in numerical computation of the inverse Mellin transform in Eq. (16) 
near z = 1. For z 1, the factor z~^ causes the integrand to vanish quickly at large \N\ 
and the integral over N converges. On the other hand, when z = 1, the convergence of the 
integral depends solely on the behavior of the Mellin-space FF pf) at large \N\. 

Since Pgg and Pgg behave asymptotically as negative constants times log \N\, while Pgq and 
Pgg vanish asymptotically as inverse powers of N, the coefficients of M_|_ and M_ in Eq. (16) 
damp the integral when as{pf) as{po). However, the integrals do not converge at z = 1 
unless fif is quite large in comparison with po- In fact, as p/ approaches po, the evolved 
FFs approximate the initial FFs, which, in some cases, are distributions at z = 1. We deal 
with this problem by rearranging the convolutions of the FFs and the PPCSs so as to treat 
the singular behavior of the FFs at z = 1 analytically. The details of the method are given 
in the Appendix. 
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FIG. 1: The ratio {da^Q/dpx)/{da-i<ii,o/dpT) for the ^Pf‘\ and channels in the process 
pp ^ H + X at ^/s = 7 TeV and \y\ < 1.2. 


IV. RESULTS FOR COMBINED LP AND NLO SHORT-DISTANCE COEFFI¬ 
CIENTS 


Now we use Eq. (5) to combine results for the LP SDCs, computed as described in Sec. Ill, 
with the SDCs through NLO in a*. For the latter, we make use of the computations in 
Refs. [1, 3, 4], taking the values of the parton distributions, rric, fir, h/, Pa, and Uf that 
are specihed at the start of Sec. III.^ 

^ In order to improve computational efficiency, we have omitted in the calculation of dcTNLO /dpx contri¬ 
butions from processes that are initiated by two light quarks, two light-antiquarks, or a light quark and 
a light antiquark, where the two initial partons can have different flavors. We use the generic expres¬ 
sion qq to denote these light-quark/antiquark initial states. The gg-initiated contributions are small in 
comparison to the sum of the qg- and gg-initiated contributions because the q and q partonic huxes are 
small in comparison to the g partonic flux. As pt increases, the sizes of the q and q partonic fluxes 
increase relative to the size of the g partonic hux because larger values of the parton momentum fractions 
are emphasized. At large values of px, dcrNLo/dpr is well approximated by da^Q/dpx- Therefore, we 
adopt the following computational strategy. In order to match what was done in the NLO calculation, 
we omit the 99 -initiated contributions in computing dcd^Q/dpx in Eq. (5). However, we take the qq- 
initiated contributions into account at large px, where they can be more important, by including them 
in the computation of da^^ jdpx in Eq. (5). Since each gg-initiated process that produces a given QQ 
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FIG. 2: The ratio {da^Q/dpr)/{da-fi-Lo/dpr) for the polarized and channels with lon¬ 
gitudinal final states in the process pp ^ H + X at y/s = 7 TeV and \y\ < 1.2. 



FIG. 3: The ratio {da^^^Q/dpx)/{da^^i^o/dpT) for the and channels in the process 

pp ^ H + X at y/s = 7 TeV and \y\ < 1.2. 
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FIG. 4: The ratio {da^^Q/dpr)/{dat^i^o/dpr) for the polarized and channels in the 
process pp ^ H + X at ^/s = 7 TeV and \y\ <1.2. h is the helicity of the QQ pair in the hnal 
state. 

We first compare our results for da^^^Q/dpT with da^i^o/dpT^ the fixed-order SDC accurate 
through NLO. Figures 1-4 show the ratios {da^Q/dpx)/{da-i<ii,o/dpr) for the polarized and 
unpolarized final states in the process pp ^ H + X at i/s = 7 TeV and \y\ < 1.2. 

In Fig. 1, we show the ratios {da^Q,/dpr)/{da-^^n^o/dpr) for unpolarized hnal states in 
the and channels. As px increases, the ratios for the and channels 

quickly approach unity because the LP-fragmentation contribution dominates the SDCs. 
This approach to unity is slower for the channel because the FF for the channel does 
not receive enhancements near z = 1 from a Dirac 6 function or plus distributions that are the 
remnants of soft divergences that cancel between real and virtual gluon-emission processes.^ 

channel contains an LP fragmentation contribution at the leading nontrivial order in cts, we can use LP 
fragmentation results to estimate the sizes of the qq-initiated contributions. These estimates indicate that 
gg-initiated contributions produce the largest fractional correction in the longitudinally polarized 
channel, in which they grow to about 5% of the total at pr = 100 GeV. Hence, we expect any errors that 
result from the omission of the gg-initiated processes in the NLO calculations to be much less than 5%. 

^ It was shown in Ref. [46] that the ratio /dpx + /dpT)/{dat^hO/dpr), which takes into account 

both the LP and NLP contributions, approaches unity much faster for the channel than does the 


13 












Pt (GeV) 

FIG. 5: The ratio /dpr)/{da-!<ji,o/dpT) for the ^Pf\ and channels with 

unpolarized final states in the process pp ^ H + X at ^/s = 7 TeV and \y\ < 1.2. 

At small Pt, the ratio {dal^^Q/dpT)/{da-!<sTo/dpT) is larger for the channel than for the 
and channels becanse, for the channel, the LO and NLO contribntions in the 
denominator tend to cancel. 

In Fig. 2, we show the ratios ((icTjjLo/'^Pr)/(<^o'NLo/<^PT) for longitndinal hnal states in 
the and channels. The approach of each of these ratios to nnity is slow. As was 
the case for the ratio of cross sections in the channel, the slow approach to nnity is a 
conseqnence of the fact that the FFs are not enhanced near z = 1 by a Dirac 6 fnnction or 
pins distribntions. 

In Fig. 3, we show the ratios {da^^^Q/dpr)/{da-!<(To/dpx) for nnpolarized hnal states in 
the color-singlet P-wave channels. We show the ratios for the polarized hnal states in 
Fig. 4. Here, h is the helicity of the QQ pair in the hnal state. The behaviors are sim¬ 
ilar to those for the ^Pf^ channel, except for the case of J = 2 with \h\ = 1, for which 
the ratio (do'NLo/'^PT)/(dcrNLo/hpT) is almost constant. We note that the deviation of 
dpT) / {daf^To/dpx) from nnity at large pt is of the same relative size as the statisti- 

ratio (dcrJi^o/^Pr)/(do-NLo/dpT). 
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FIG. 6: The ratio /dpT)/{da^i,o/dpT) for the polarized and channels with 

longitudinal final states in the process pp ^ H + X at ^/s = 7 TeV and |y| < 1.2. 


cal uncertainty in the NLO calculation. 

Next we compare /dpr, the SDC that includes both the hxed-order cor¬ 

rections through NLO and the additional LP corrections, with duNLo/dpT, the SDC 
that includes hxed-order corrections through NLO. Specihcally, we show the ratios 
(dcr^Lo^^*^/^PT)/(do'NLo/dpT) for the polarized and unpolarized hnal states in Figs. 5-8. 
With the exception of the channel, the additional LP-fragmentation contributions are 
of the order of 100% at large pt- As we have mentioned, because there is a partial cancel¬ 
lation between the LO and the NLO contributions in the channel, the additional LP 
fragmentation corrections have a signihcant impact on the shape in that channel. For the 
channel the additional LP-fragmentation contributions are negative and only mildly 
alter the shape. 

As was pointed out in Ref. [25], the effects from the all-orders resummation of logarithms 
of Pr/ml are small. In the case of the channel, almost all of the effects of the large 
logarithms are already accounted for in the NLO contribution. In the cases of the and 
channels, the all-orders resummations of logarithms shift the FFs by only about 2% 
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FIG. 7: The ratio /dpT)/{da^i,o/dpx) for the and channels with unpolarized 

final states in the process pp ^ H + X at -y/i = 7 TeV and \y\ < 1.2. 



PT (GeV) 


FIG. 8: The ratio /dpx)/{daj^i^o/dpx) for the and channels with polarized 

final states in the process pp ^ H + X at y/s = 7 TeV and \y\ < 1.2. 
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and 5%, respectively, at pt = 52.7 GeV because contributions from the running of as and 
from the DGLAP splitting cancel. Hence, almost all of the large additional LP corrections 
that we hnd arise from nonlogarithmic contributions of order a^. 

Finally, we discuss the LP-fragmentation contribution to the channel. Since the 
FF for this channel begins at order the channel receives an LP contribution that 
begins at order (NNLO). We do not include this LP-fragmentation contribution in our 
analysis. However, we have estimated its size by making use of the FF at order a^. At 
Pt = 10 GeV, the LP contribution is about an order of magnitude smaller than the hxed- 
order contribution through NLO. The LP contribution reaches the same size as the hxed- 
order contribution through NLO at around px = 50 GeV. Finally, when px = 130 GeV, the 
LP contribution is almost an order of magnitude larger than the hxed-order contribution 
through NLO. Although the LP-fragmentation contribution can have a signihcant effect on 
the color-singlet contribution at large px, its effect on the cross section is only of the order 
of 1% of the measured cross section at px = 130 GeV. 

V. FITS OF CROSS-SECTION PREDICTIONS TO DATA 

In this section we extract the color-octet LDMEs by htting the cross-section predictions 
that are based on the LP-I-NLO SDGs to the measured cross sections. We use the resulting 
LDMEs to make predictions for the prompt-J/'^ polarization. In order to suppress possible 
nonfactorizing contributions, we ht only to data for which px is greater than 3m//, where 
thh is the quarkonium mass. Since the shape of the px distribution determines the LDMEs, 
it is crucial to use the data at the highest pr values in the hts. 

In the case of the direct J/ijj cross section, we estimate the theoretical uncertainties in 
the SDGs to be 25% of the central values. We arrived at these uncertainties by varying 
the factorization scale pf and the renormalization scale pr independently between \mx and 
2mx- This 25% uncertainty is also roughly the size of the uncertainty that one would expect 
from uncalculated corrections of higher order in v. In the cases of the cross sections of the 
excited charmonium states, we take the uncertainties to be 30% of the central values because 
the for those states is larger than for the J/'^. 
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A. Production of 'ip(2S) 


We determine the three color-octet LDMEs by performing a least-y^ fit to the 

CDF [47] and CMS [12, 48] cross-section data. In order to suppress possible nonfactorizing 
contributions, we use only the data for which is greater than 11 GeV. We ignore feeddown 
contributions from decays of heavier quarkonia. 

In the case of the color-singlet LDME, we take a value that was determined in a potential- 
model calculation [49]: = 0.76 GeV^. Different choices for the value of the 

color-singlet LDME would have little effect on our results, as the contribution from the 
color-singlet channel is much smaller than the theoretical uncertainties. In the lowest pt 
bin that we consider for the CMS data that have \y\ < 1.2 (11 GeV < < 12 GeV), the 

contribution from the color-singlet channel is only about 5% of the cross section, and the 
color-singlet contribution drops to 0.2% in the highest pt bin (75 GeV < pt < 100 GeV). 

The htted LP-I-NLO cross section is compared with the data in Fig. 9. The quality of the 
ht is quite good, with y^/d.o.f. = 1.71/29. As can be seen from Fig. 10, the cross section is 
dominated by the channel at moderate values of but not at large values of pt- The 
concept of dominance has been suggested previously in Refs. [4, 25, 50]. 

The color-octet LDMEs that are obtained from the ht are 


(C>’^(2^)(34®1)) = (-1.57 ± 2.80) X 10-3 GeV3, (17a) 

{C>’^(2^)(i5'fl)) = (+3.14 ± 0.79) X 10-2 GeV3, (17b) 

(C>^( )( Pq )) ^ 4. ^ 21) X 10-3 GeV3. (17c) 

mi 


The uncertainties that are shown above are correlated. The correlation matrix of the 
uncertainties in {C>^*^^'^)(3S'f^)), ((^^(^'^^(^S'o®')), and Pf^))/ respectively, is 


^ 7.85 -14.7 3.36^ 


(jp{2S) _ 


\ 


-14.7 62.2 -5.52 

3.36 -5.52 1.46 


X 10-^^ GeV®. 


(18) 


It is useful to examine the correlation matrix of relative uncertainties, whose 

components are dehned by 


(jp{2S) _ 


^i>{2S) 

^nm 

o o 

r)\^ m. 


(19) 
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FIG. 9: The differential cross sections for prompt tpi^S) production at the Tevatron {^/s = 
1.96 TeV) and the LHC {^/s = 7 TeV). H^( 25 ) = Br[V'(2S') —^ where Br denotes a branching 

ratio. 


where On is the central value of the nth LDME. Then is given by 


(f7h(25) 


^32.0 2.98 18.9^ 
2.98 0.63 1.55 
^18.9 1.55 11.3^ 


X 10"^ 
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FIG. 10: Contributions of the individual channels to the prompt ^^{23) differential cross section 
at the LHC (v^ = 7 TeV). = Br[V’(25) ^ 


The normalized eigenvectors of are 

^ 0.858 ^ 


M2S) 


0.0780 

0.508 


/ n 1 Qi \ 


M2S) 


\ 

M2S) 


0.181 

0.879 

-0.441 


/ O /IQ1 \ 


M2S) 




\ 


-0.481 

0.470 

0.740 
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/ 


and the corresponding eigenvalues are = 4.34, = 4.67 x 10“^, and = 

1.96 X 10“^. The eigenvector is predominantly and its uncertainty [(A 2 

is fairly small. On the other hand, the eigenvector has a very large uncertainty 

[(Ay( 2 S))i/ 2 ]^ Hence, the and LDMEs can vary together in a correlated way that 
tends to preserve the dominance. (Recall that the SDCs for these channels have opposite 
signs.) The eigenvector has a very small uncertainty and, therefore, the 

anticorrelated variation of the and ^Pq*^ LDMEs is highly constrained. 


ih(25) 


B. Production of Xci and Xc2 

We determine the two LDMEs for Xcj by htting to ATLAS cross-section data [51]. In 
order to suppress possible nonfactorizing contributions we £t only to data for which pt is 
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FIG. 11: The differential cross sections for prompt Xci and Xc 2 production at the LHC (+5 
7 TeV). = Br[xcj ^ J/^ + 7 ] x Br[J/'!/’ -)■ 


greater than 11 GeV. We ignore feeddown contributions. The 'ip{2S) decays into Xcil and 
Xc 27 with branching ratios of 9.55% and 9.11%, respectively. These contributions amount to 
only a few percent of the measured cross sections and are much smaller than the theoretical 
uncertainties. 

The htted LP+NLO Xci and Xc 2 cross sections are compared with the data in Fig. 11. 
We do not consider the Xco cross section because the Xco branching ratio to J/^jx is small 
and the corresponding contribution to the prompt J/ijj cross section is negligible. Again, 
we obtain a good ht to data, with y^/d.o.f. = 1.19/8. The contributions of the individual 
channels to the prompt-Xci and prompt-Xc 2 cross sections are shown in Fig. 12. There are 
substantial cancellations between the contributions of the and channels. 
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FIG. 12; Contributions of the individual channels to the differential cross sections for prompt Xci 
and Xc 2 production at the LHC (-v/i = 7 TeV). = Br[xcj —>■ J + 7 ] x Br[J/V' —)• 


The resulting LDMEs are 


((^xc(3^[8])) _ (5.74 ± 1 . 31 ) X 10-3 GeV3, 
^ = (3.53 ± 1.08) X 10-2 GeV3. 




(22a) 

(22b) 


The correlation matrix of the uncertainties in ((^^‘'(^S'f^)) and respectively, 


IS 


, 1.71 14.0 , 

= \ lx 10-® GeV®. 

14.0 117 


(23) 


The relative uncertainties in these LDMEs are fairly small, but there are substantial corre¬ 
lations between them. The correlation matrix of relative uncertainties is 


= 


5.18 6.91 
6.91 9.39 


X 10 


-2 


(24) 


The normalized eigenvectors of are 


vf = 


0.595 

0.804, 




0.804 

-0.595 


(25) 
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and the corresponding eigenvalues are = 0.145 and A^'" = 5.70 x 10“^. We see that, 
while the eigenvector Vi'" has a small uncertainty, the and LDMEs can vary in a 
correlated way. However, from the very small uncertainty of the eigenvector we see that 
anticorrelated variation of these LDMEs is highly constrained. 

At leading order in v, the color-singlet LDME is related to the derivative of the wave 
function at the origin as 

=2Af,4^|fl'(0)|^ (26) 

The value of the color-singlet LDME that we obtained from our £t corresponds to |i?'(0)p = 
0.055 ± 0.017 GeV®. This is consistent with the value |i?'(0)p = 0.075 GeV^ that was 
obtained in Ref. [49] by using the Buchmiiller-Tye potential. It is also consistent with the 
value that was determined in Ref. [52] from the two-photon decay rates of the Xco and the Xc 2 , 
namely, = 0 . 060 lo;o 29 GeV^, which corresponds to |R'( 0 )p = 0.042to;o2o GeV®. 


C. Production of prompt J/V’ 


We determine the J/ip LDMEs by fitting to the GDF [11] and GMS [12, 48] prompt-J/'^ 
cross-section data. In order to suppress possible nonfactorizing contributions, we fit only to 
data for px greater than 10 GeV. We compute the feeddown contributions from the decays 
of ip{2S), Xci, and Xc 2 by making use of the LDMEs that were determined in the preceding 
sections. The prompt-J/ ip cross section is given by 


dpT 


d(^j/^ + ^^^^Br[^(2A) ^ J/iP + X] + ^Br[x,i ^ J/iP + 7 ] 


dpT 

da 




dp. 


dp: 


y^Br[xc2 —t J/ip 


(27) 


Here, we ignore the feeddown contribution from the decay of the Xco- As we have mentioned, 
the XcQ decays into with a small branching ratio, and the contribution to the prompt 

J/ip cross section is negligible, pt is the transverse momentum of the J/ip^ and is the 

transverse momentum of 77 = ip{2S), Xcj- In the feeddown contributions, we take to be 


H 

Pt = - Pt- 


(28) 


23 



The relation (28) is derived by neglecting the 3-nionientuni of the J/'ip in the H rest frame 
in comparison with 

We take the value of the color-singlet LDME that has been obtained from the electro¬ 
magnetic decay rate [55]: = 1.32 GeV^. Again, the contribution from the 

color-singlet channel is much smaller than the theoretical uncertainties, ranging from 4% for 
the bin 10 GeV < pr < H GeV to 0.2% for the bin 95 GeV < pr < 120 GeV in comparison 
with the direct J/ip cross section. 

We obtain a good £t to the data, with y^/d.o.f. = 8.20/40. The fitted LP-I-NLO cross 
section is shown in comparison with the data in Fig. 13. The contributions of the individual 
channels to the direct J/ip cross section are shown in Fig. 14. The direct J/ip cross section 
is dominated by the channel at all values of pt between 10 GeV and 100 GeV. 

The color-octet LDMFs that are obtained from the fit are 


((^j/V>(3^[8])) _ ^ ^ 3Q-3 Qgys^ ^29a) 

^ QJ4) X ^Q-i GeV^ (29b) 

(C)J/b(3p^ )) ^ ^ ^ 

mj 


The correlation matrix of the uncertainties in and 

((!l'^/'^(3P(|®]))/m^, respectively, is 


(jJ/tp 


^ 13.3 -38.2 5.48^ 
-38.2 188 -14.6 

^ 5.48 -14.6 2.29^ 


X 10"^^ GeV®. 


(30) 


^ We have estimated the effects of corrections to this relation on the contributions of Xci and Xc 2 feeddown to 
the J/b unpolarized and polarized cross sections. In these estimates, we computed the angular distribution 
of the J/ip momentum in the Xcj rest frame by making use of the formalism of Ref. [53], and we included 
the El, M2, and E3 electromagnetic transition amplitudes, taking the M2 and E3 amplitudes to be 
given by the central values of the measurement of the CLEO Collaboration [54]. We find that the 
corrections to the feeddown contributions are no more than 8% in any of the J/ip polarization channels. 
Furthermore, the corrections are essentially flat as functions of pt, deviating by only about 1% over the 
range 10 GeV < pr < 100 GeV, with almost all of the deviation occurring between 10 GeV and 15 GeV. 
Hence, the corrections have little effect on the shapes of the cross sections and can be absorbed into 
normalization shifts of the LDMEs of a few percent or less. 
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FIG. 13: The differential cross section for prompt J/'ip production at the Tevatron {y/s = 1.96 TeV) 
and the LHC {y/s = 7 TeV). Bjj^ = Br[J/V’ —)• 


The correlation matrix of relative uncertainties is 




^26.1 4.88 24.7^ 
4.88 1.55 4.26 
^24.7 4.26 23.5^ 


X 10"^. 


(31) 
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FIG. 14: Contributions of the individual channels to the differential cross section for direct J/V’ 
production at the LHC (-y/s = 7 TeV). = Br[J/'0 —)• 


The normalized eigenvectors of are 


j/i/> 

Vi 


( 0.719^ 
0.131 
yO.682 j 


^ 0.168^ 


= 


0.920 

-0.354 


^^- 0 . 674 ^ 


= 


V 


0.369 

0.640 


(32) 


and the corresponding eigenvalues are = 0.504, = 8.06 x 10“^, and = 

2.35 X 10“^. As is the case for the ■^(25'), the eigenvector that is predominantly 
namely, has a fairly small uncertainty. However, the eigenvector has a very large 
uncertainty. Therefore, variations of the and LDMEs are correlated and tend to 

preserve the dominance. (Recall that the SDCs for these channels have opposite signs.) 
The very small uncertainty of the eigenvector means that the anticorrelated variation 
of the and LDMEs is highly constrained. 
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FIG. 15: Fraction of prompt J/ijj's produced in feeddown from ^p{2S) decays at the LHC {^/s 
7 TeV). 
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FIG. 16: Fraction of prompt J/tp’s produced in feeddown from Xci and Xc 2 decays at the LHC 
(v^ = 7TeV). 
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VI. PREDICTIONS FROM EXTRACTED LDMEs 


In this section, we use the LDMEs that we have extracted from the hts to cross sections to 
make predictions of cross-section ratios and polarizations. We estimate the uncertainties in 
these predictions by making use of the eigenvectors and eigenvalues of the LDME uncertainty 
correlation matrices. In the expression for each prediction, we write the LDMEs in terms of 
the eigenvectors. Then, we vary each eigenvector about its central value by an amount that 
is equal to the square root of its eigenvalue. We take the resulting variation in the prediction 
as the uncertainty in the prediction from variations of that eigenvector. Finally, we estimate 
the total uncertainty in the prediction by adding the uncertainties from the variations of the 
individual eigenvectors in quadrature. 


A. Ratios Rh 


We can use our predictions for the J/t/’, 'ip{2S), and Xcj cross sections and the LDMEs 
that we have extracted to compute the ratios Rh, which are dehned by 


Rh 


Br[H J/ij + X] X dan/dpif 


da™* 


/dpT 


(33) 


where pif is given in Eq. (28). In Figs. 15 and 16 we show our results for Rp( 2 S) and 
R^^ = R^^^ + Rxc 2 ) respectively. As can be seen from Fig. 16, our prediction for R^^ lies 
systematically below the ATLAS [51] and LHCb [56] measurements for pt < 15 GeV. This 
discrepancy occurs because the prediction for the numerator of R^^ lies slightly below the 
data at low pr, while the prediction for the denominator of R^^ lies slightly above the data 
at low Pt- However, the predictions for both the numerator and the denominator agree with 
the data within uncertainties. We also note that corrections to the relation (28) for the .J/ijj 
momentum would increase the theoretical prediction for R^^ by a few percent. 


B. Polarization predictions 

We now compute prompt-'^(2S') and prompt-J/'^ polarizations by making use of the 
LDMEs that we have determined from hts to the cross-section data. We also compute the 
effects of feeddown from the t/’(25') and Xcj states on the polarizations of the prompt J/'ip^s. 



For J = 1 states, one measure of the polarization is the polarization parameter Ag, which 


is dehned as 


a — 3ar 


-^e - -^ 

+ a L 

where a and ai are the polarization-summed and longitudinal cross sections, respectively. If 
the J = 1 state is completely transversely (longitudinally) polarized, then az, = 0 (a/, = a), 
and Xg = -|-1 {Xg = —1). If the J = 1 state is unpolarized, then a = Saz,, and Xg = 0. 

We show the polarization of the as produced at the LHC at i/i = 7 TeV and at 

the Tevatron at y/s = 1.96 TeV in Figs. 17 and 18, respectively. The prediction for the CMS 
polarization is in fair agreement with the CMS data [15]. The prediction for the polarization 
at the Tevatron is in rough agreement with the CDF Run I [13] and Run II [14] data, given 
the very large error bars. (Although the CDF Run I data were taken at i/i =1.8 TeV, rather 
than at i/i = 1.96 TeV, this energy shift produces a negligible change in the polarization 
prediction.) The predicted t/’(25') polarization grows as px increases, owing to the fact 
that the channel is no longer dominant at large px- Hence, measurements of the 
polarization at larger values of px would provide an important test of the theoretical 
prediction. 

The longitudinal prompt-J/'^ cross section, including the feeddown contributions from 
the decays of the 'ip{2S), the Xci, and the Xc 2 is computed as follows: 


'^'^J/y(A=0) _ '^'^J/y(A=0) d(T^(25)(A=0) 

dpx dpx dpp-‘^^^ 


Br[^(2A) ^ J/ij + X] 




2 d(7^^2(A=0) ’^*^Xc 2(A=+1) *^‘^Xe2(A=-l) 

3 dpf^ 2 V dpf^ dpf^ 


Br[xc 2 + 7], (35) 


where Px is given by Eq. (28). In deriving Eq. (35), we have assumed that the polarization 
of the t/’(25') is completely transferred to J /ip and that the decays Xcj —t J / ip + x proceed 
through an El transition. In the decays Xcj J/'P't'Xi fhe higher multipole corrections are 
poorly known, but they have little effect on the polarizations of the J /ip’’s that are produced 
in Xcj decays [53]. 

We show the polarization of J/ip^s from Xcj decays at the LHC at = 7 TeV in 
Fig. 19. In Fig. 20, we show the polarization of prompt J/ip^s produced at the LHC at 
a/s = 7 TeV, including feeddown from the ip{2S) and the Xcj states. The prediction is 
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FIG. 17: Polarization of prompt '0(25') at the LHC = 7 TeV). 
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FIG. 18; Polarization of prompt 0(25) at the Tevatron (-^/i = 1.96 TeV). 
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FIG. 19: Polarization of J/'ijj from Xcj decays at the LHC {^/s = 7 TeV). 

in good agreement with the CMS data [15]. Finally, in Fig. 21, we show the polarization 
of prompt J/'ip^s produced at the Tevatron at y/s = 1.96 TeV, including feeddown from 
the 'ip{2S) and the Xcj states. The prediction is in good agreement with the CDF Run I 
data [13], but disagrees with the CDF Run II data [14]. (Although the CDF Run I data 
were taken at \/s = 1.8 TeV, rather than at ^/s = 1.96 TeV, this energy shift produces a 
negligible change in the polarization prediction.) We note that the predicted polarizations 
are almost the same for the LHC and the Tevatron, while the CDF Run II polarization data 
lies signihcantly below the CMS polarization data. 

The fairly small polarizations that are seen in the predictions for the prompt J/'ijj^s and 
'^(25')’s are a consequence of the dominance in the production rates of the channel, 
which, of course, is completely unpolarized. This mechanism whereby small polarizations 
can be obtained was noted previously in Refs. [4, 25, 50]. 
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FIG. 20: Polarization of prompt at the LHC {^/s = 7 TeV). The polarizations of 

produced in feeddown from the 'tp{2S) and XcJ states are shown with dashed and dotted lines, 
respectively. 
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VII. SUMMARY AND DISCUSSION 


In this paper, we have computed, in the NRQCD factorization framework, leading-power 
(LP) fragmentation corrections to production of the charmonium states J/ip, Xcj, and 
in pp collisions at the Tevatron and in pp collisions at the LHC. Specihcally, our calculation 
makes use of parton production cross sections (PPCSs) through order af (NLO) and frag¬ 
mentation functions (FFs) through order We have also used the DGLAP equation to 
resum leading logarithms of p^/m?^ to all orders in Our calculations take into account 
the effects of feeddown from the 'ip{2S) and Xcj states on the prompt-J/V’ cross sections 
and polarizations. Hence, the work in the present paper is an extension and a rehnement 
of the work in Ref. [25], which also addressed LP corrections, but which did not include 
computations of cross sections or polarizations for the '0(2*5') or XcJ states or include the 
effects of feeddown from those states. We hnd that the LP corrections, beyond those that are 
contained in hxed-order calculations through NLO in as, are substantial—typically of order 
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FIG. 21: Polarization of prompt J/V^’s produced at the Tevatron {-^/s = 1.96 TeV). The polariza¬ 
tions of J/V^’s produced in feeddown from the V'(2S') and Xcj states are shown with dashed and 
dotted lines, respectively. 

100% at large px- Owing to a partial cancellation between the LO and NLO contributions 
in the channel, the LP corrections have a very signihcant effect on the shape in that 
channel. 

As was pointed out in Ref. [25], the all-orders resummations of logarithms of 
have only small effects on the predictions for the cross sections and polarizations. Hence, 
almost all of the large additional LP corrections that we hnd arise from nonlogarithmic 
contributions of order a®. 

Our approach in calculating the LP fragmentation corrections is to use the most accurate 
results for the PPCSs and FFs that are currently available. This means that we have 
computed some, but not all, of the LP contributions in order af, i.e., NNLO in terms of 
the hxed-order calculations. In the case of gluon fragmentation, a complete calculation of 
the order-a® contributions in the and channels would require the calculation of 
the NLO corrections to the FFs for those channels. In the case of gluon fragmentation, 
a complete calculation of the order-af contributions in the channel would require a 
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calculation of the NNLO corrections to the FFs for that channel and the NNLO corrections 
to the PPCSs. We expect these uncalculated LP corrections in order to be of comparable 
size, in each channel, to the LP corrections that we have calculated in this paper. Hence, 
the actual theoretical uncertainties may be much larger than the estimates that we have 
obtained by varying the scales and fif. However, we emphasize that the calculation in 
this paper eliminates the largest existing source of theoretical uncertainty by taking into 
account leading logarithms of at all orders in o;*. 

We have combined the LP corrections that we have calculated with the NLO fixed-order 
calculations from Refs. [1, 3, 4] to obtain predictions for the production cross sections and 
polarizations as functions of pt- By htting the cross-section predictions to the Tevatron 
and LHC cross-section data, we have obtained values for the NRQCD nonperturbative long¬ 
distance matrix elements (LDMEs) that enter into the production predictions through order 
Since the LP approximation is valid only for pt 3> mn, where mn is the quarkonium 
mass, we use only data for which pr is greater than SniH- We obtain good hts to the high-p^ 
cross sections, with y^/d.o.f. 1 in each case. 

One interesting result of the hts to the Xcj cross sections is that the value of the 
LDME that we obtain is in good agreement with the value that has been obtained in a 
potential model and with values that have been extracted from the two-photon decays 
of the Xco and Xc 2 - This agreement of values of the ^Pq^^ LDME that have been obtained 
through very different methods is important evidence in support of the NRQCD factorization 
conjecture. We note that, in previous works on the Xcj cross section, which were based on 
hxed-order NLO calculations, the ^P^^^ LDME was hxed to values that were obtained from 
potential models [6, 57-61]. 

We have used our cross-section predictions to predict the ratio which is the Xcj 
feeddown contribution to the prompt-J/'^ cross section divided by the prompt-J/'^ cross 
section itself. The prediction lies systematically below the data for px < 15 GeV. This 
discrepancy in seems to be the result of a downward deviation in the numerator combined 
with an upward deviation in the denominator. However, the predictions for both numerator 
and the denominator agree with the data within uncertainties. 

We have also used the extracted LDMEs to predict the J/t/’, 'ip{2S), and Xcj polarizations. 
The predictions for the J/'ip polarizations agree with the CMS data and the CDF Run I 
data, but he systematically above the CDF Run H data. The CDF Run I data show a 


34 



slightly longitudinal polarization, while the CMS data show a slightly transverse polarization. 
However, the theoretical predictions are very similar for the CDF and CMS kinematics. 
The predictions for the polarizations agree with the Tevatron data and the LHC 

data, although the theoretical and experimental uncertainties are quite large. For 'ip{2S) 
production, the channel is no longer dominant at large pt, and, so, the predicted 
polarization becomes more transverse as pt increases. It is important to test this 
prediction through measurements of the polarization with good precision at larger 

values of pt- There are, as yet, no measurements of the XcJ polarizations. These would also 
provide very useful tests of the theoretical predictions. 

While we have obtained a reasonably good description of the hadropro duct ion cross sec¬ 
tions and polarizations for the J/^j, 'ip{2S), and Xcj states, our results do not address two 
outstanding problems in quarkonium production, namely, the HERA J/ijj photoproduction 
cross section, as measured by the HI Collaboration [16, 17], and the pc hadroproduction 
cross section, as measured by the LHCb Collaboration [20]. In the case of the J/ijj photo¬ 
production cross section, additional LP fragmentation corrections, analogous to those that 
were computed in this paper, were computed in Ref. [62]. Those additional LP corrections 
have very small effects on the photoproduction cross section. For the choices of LDMEs 
that were used in Ref. [62], the theoretical prediction for the photoproduction cross section 
is dominated by the contribution from the channel. The value for in 

Eq. (29) is about 10% larger than the value from Ref. [25], which was used in Ref. [62]. 
Hence, it makes the discrepancy between theory and experiment slightly worse. In the case 
of the Pc cross section, the change in value of from Ref. [25] to the present 

paper also makes the discrepancy between theory and experiment slightly worse. 

While there remain important discrepancies between theory and experiment in quarko¬ 
nium production at high pt-, the theoretical predictions are far from settled. At a minimum, 
a complete calculation of all of the LP contributions in order af is needed in order to have 
reasonable control of the theoretical uncertainties. These LP contributions in order may 
be most important in the channel because of their greater potential to affect the shape 
in that channel. Higher-order calculations of NLP contributions may also be needed, es¬ 
pecially in the channel, for which the LP contributions are not dominant until very 
large values of pt- New measurements at the LHC of the 'ip{nS), Xcj, ^{nS), and Xbj cross 
sections and polarizations and the pc cross section, all at unprecedentedly large values of pt, 
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can provide definitive tests of the improved theoretical predictions. 
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Appendix: Numerical treatment of divergences in fragmentation functions 


The LP-factorization contribution to the cross section is given by the convolution of the 
PPGSs daAB^i+x/dpT and the FFs pf): 


AB^QQ{n)+X 

dpT 


= [ dz 
J ZQ 


da 


AB^i+X 




dpT 


A 


i^QQ(n) \ 


b h/)- 


(A.l) 


Here, Zq = ^(e+^ + e“^). We compute the evolved FFs by solving the LO DGLAP equation 
in Mellin space (moment space) and performing the inverse Mellin transform numerically. 
As is discussed in Sec. HI G, the inverse Mellin transform becomes numerically unstable 
near z = 1 because pf) can vary rapidly in this region and may even diverge at 

z = 1. 
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In order to deal with this problem, we partition the integral over z as follows: 
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where e is a small, positive number that is chosen so that the evolved FF Di^QQ(ri){.Z) Pf) 
can be computed reliably for < 1 — e. In order to compute the integral over 1 — e < < 1, 

we use the fact that the PPCSs behave as daAB^i+x/dpr ~ for z ^ 1, where N ^ A. 
Hence, we have 
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, (A.3) 


where we have expanded z~^ddAB^i+x/dpr in powers of 1 —x and retained only the leading- 
order contribution, which is simply the value of ddAB^i+x/dpr at z = 1. The quantity 
+ 1)/^/)) {di + l)st moment of is known analytically, and the integral 

over the range 0 < x < 1 — e can be computed numerically. Hence, in our calculations, we 
use the expression 
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(A.4) 


In our numerical calculations, we take iV = 4 and e = 10“®. We have varied N between 3.1 
and 7 and hud that the largest sensitivity to N occurs at low pt and is less than 3 x 10“® 
of the contribution in each channel. 
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We have compared numerical results from Eq. (A.4) for pj near po with the analytic 
expression for the evolved FFs through second order in a.^- The results agree to better than 
1%. We expect numerical difficulties in Eq. (A.4) to be most severe as p/ approaches po; 
where the evolved FFs approximate the initial FFs, which, in some cases, are distributions 
at z = 1. Hence, good agreement with the analytic expressions in this region gives us 
conhdence that the algorithm that is based on Eq. (A.4) is reliable. 
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